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ABSTRACT 



The underlying purpose of this research was to understand 
and predict the response of upper ocean boundary layer 
mixing to high frequency oscillations in the atmospheric 
forcing attributable to the diurnal heat flux cycle, 
unsteadiness in the wind, and other short term changes. To 
accomplish this task, a non-stationary , one-dimensional 
bulk model of the mixed layer originally proposed by 
Garwood (1976, 1977) is generalized by deriving a new 
equation for the entrainment buoyancy flux which includes 
the unsteady term. To examine the importance of the un- 
steady term, a reduced form of the turbulent kinetic energy 
equation is solved numerically. The results predict the 
high frequency cutoff above which the unsteady term should 
not be neglected. The quasi-steady state assumption is 
accurate only for low frequency forcing up to this cutoff. 
There was also some dependence upon the amplitude of the 
imposed forcing. To investigate the asymptotic properties 
of the mixed layer model, the full form of the Garwood 
model is solved numerically and compared with the solutions 
to the de Szoeke and Rhines model. The results of this 
research provide a theoretical basis for realistically 
applying the mixed layer models to ocean boundary layer 
problems on all time scales. 
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NOTATION 



) - lim ± f 

^ ’ - U<o >-* hy^ As/x 


( ) dy Horizontal mean 


<( IT J-h ^ ^ 


Vertical mean over the 
mixed layer 


a 


Amplitude of fluctuation 


b = b + b’ 


Buoyancy 


C = C + c' , CxU +!'<' 


Velocity in complex form 


Q 


Specific heat at con- 
stant pressure 


d 


(1) Ekman depth of 
frictional resistance 
for the upper ocean 

(2) Depth of surface 
energy perturbation 
zone 


E = a'^+ V'* + w'^ 


Total turbulent kinetic 
energy per unit mass 


f = 2 ini SiTJ 4 


Coriolis parameter 


g 


Apparent gravitational 
acceleration 


h 

-l4'vV7o)h 

” — E — 


Depth of mixed layer 

The mixed layer 
stability parameter 


L 


Obukhov length scale 


1 


Size of the large eddies 


m, m. 


Constants 


N 


Brunt-Vaisala frequency 


P = pg + P' 


Pressure 
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Geostrophic pressure 
component 


Q 


Radiation absorption 


Ri 


Bulk gradient Richardson 
number 


az/ az* 


Local gradient Richardson 
number 


^ LI# 

T = T + T' 


Frictional Rossby number 
Temperature 


Tp 


Nondimensionalized 
forcing period 


t 


time 


Uj_ = (u +u' 1 V + v{ w +wM 


Total instantaneous 
velocity (without 
geostrophic component) 


u* 


Frictional velocity 


Ug« 


Geostrophic wind velocity 


w 


Wind velocity 


Xj_ = ( Xi y. z) 


Rectangular space 
coordinates 


oL=:±^ 

P* 9T 

r = N* = -^ for 2 <-h-<f 

92. 


Thermal expansion 
coefficient 


T 


Lapse rate of potential 
temperature above the 
inversion base 




Thickness of entrainment 
zone 


r Co 


Excess surface velocity 


A'&s b(-h) -b (-h-<T) 


Change in mean buoyancy 
across the entrainment 
zone 
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T 

1 T(t)i =ie.oi 
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rp 
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I. INTRODUCTION 



A. PURPOSE OF THE STUDY 

The objective of this research was to study the one- 
dimensional mixed layer model of Garwood (19?6, 197?) in 
terms of a possible generalization of the extrainment 
buoyancy flux equation and to examine the responses of the 
mixed layer turbulent energy intensity to the different 
amplitude and frequencies for wind-shear turbulence prod- 
uction. Because of the frequency dependence, the steady 
state ( , where <li> is the turbulent kinetic energy 

averaged over the mixed layer depth) or quasi-steady state 
( —4- ^" 0 ) assumption is good in some frequency ranges 

but it is not possible in other ranges. Using a constant 
dissipation time scale, Garwood (19?6) suggested that if 
the quasi-steady state assumption is made to facilitate the 
solution without filtering the surface boundary conditions, 
an incorrect high frequency response will not only be pres- 
ent but may bias the mean mixing trend. 

A further objective of this study is to critically 
examine the asymptotic time regimes of the Garwood m.odel 
in comparison with those of the de Szoeke and Rhines (1976) 
model. These two tasks will provide a theoretical basis for 
realistically applying the mixed layer model to ocean 
boundary layer problems on all time scales. 



• y 



B . METHOD 



The equations of the Garwood (1976, 1977) model are 
listed and their derivation is discussed briefly. The 
de Szoeke and Rhines (1976) work is also reviewed in 
order to provide the background information necessary to 
compare the asymptotic properties of the two different 
models. To generalize the Garwood entrainment buoyancy 
flux equation, the unsteady term of the turbulent kinetic 
energy equation at the base of the mixed layer, which was 
usually considered to be negligible by earlier investigators 
is parameterized. For the study of frequency-dependent 
response to atmospheric forcing, the method applied by 
Garwood (1976) in which he used a reduced and nondimens ional 
ized form of the turbulent kinetic energy equation, and 
solved analytically by using a constant dissipation time 
scale, is reviewed. This same form of the equation is 
solved numerically by using a non-constant dissipation time 
scale, » provided by Garwood (1976), where h is the 

mixed layer depth. In this case, the solutions for non- 
dimens ional ized frequencies 2oott , , .2tt , 

and .02TT , where u; is the angular frequency of the 

normalized wind forcing and u^ is the maximum friction 
velocity, and amplitudes of fluctuations a = 0.1. 0.5 and 

1 for the wind production are compared. 

To study the asymptotic properties, the Garwood model 
is first solved numerically for different stabilities and 
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its results are compared with those of de Szoeke and Rhines 
by using equivalent calibration constants, boundary con- 
ditions and initial conditions. 

C. HISTORICAL LITERATURE REVIEW 

The vertical fluxes of heat, salt and momentum across 
the sea-air interface are the source of almost all oceanic 
motions. In the fully turbulent oceanic mixed layer 
bounded by the sea-air interface above and by the dynamically 
stable water mass below, the vertical fluxes are large. 

Below the mixed layer they are usually negligibly small so 
that one can decouple the mixed layer from the underlying 
stable, quiescent water mass. The immediate local reaction 
of this mixed layer to those fluxes result in a homogeneous 
water column, i.e., vertical uniformity in the mean velocity 
and density. This adjustment of the density and velocity 
structure of the surface layers of the ocean to variable 
fluxes has been the subject of a large number of studies 
since Ekman (1905) 's treatise where he originated the 
concept of a depth, d = 7 • 5csec'*) is wind velocity 

and 4 is latitude), of frictional resistance for the upper 
section of a wind stressed ocean and suggested that any 
surface mixed layer by the action of wind has the same 
order of this depth. This depth comes from the mathematical 
solution to the steady state horizontal momentum equation 
in which Reynold's stress is related to the mean shear by 
a constant eddy viscosity, K . Rossby and Montgomery ( 1935 ) 
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pointed out that the depth of a surface drift current layer, 
h ccU*/-f ( I is Coriolis parameter) , and Ekman depth are not 
necessarily comparable. Both depths are derived by con- 
sidering only the momentum budget, neglecting the effects 
of buoyancy and the mechanical energy budget. 

Munk and Anderson ( 19 ^ 8 ) developed a simultaneous 
solution of steady-state heat and momentum profiles by 
using a viscosity and an eddy conductivity that depend upon 
the gradient Richardson number. They did not recognize the 
presence of a sharp interface at the bottom of the mixed 
layer and therefore their result is closer to Ekman 's 
solution than to the physical reality. Kitaigorodsky (i960) 
concluded by dimensional analysis that the mixed layer depth 
must be proportional to the Obukhov length scale, L. The 
basic flaw in his steady state model is the assumption that 
the entire ocean mixed layer is analogous to the constant- 
flux atmospheric surface layer, because with his constant 
flux at the surface the mixed layer temperature and depth 
cannot both remain unchanged. 

Kraus and Rooth (I96I) developed a plausible steady state 
model based primarily upon the buoyancy equation. In their 
model, steady state is achieved by balancing solar radiation 
with a net surface heat loss due to evaporation, conduction 
and long wave back radiation. The unstable density profile 
above the compensation depth is a convective source of the 
turbulent kinetic energy for mixing. Recognizing the limi- 
tation in Kraus and Rooth model--no provision for a possible 
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dovmward surface heat flux, no account of mechanical 
production of turbulent kinetic energy and the steady state 
constraint — Kraus and Turner (I 967 ) improved the one- 
dimensional model hy considering the turbulent kinetic 
energy budget utilizing the two separate equations, the 
heat equation and the mechanical energy equation. Because 
the frictional generation of heat is negligible, the 
vertically integrated heat equation provides a relationship 
for the conservation of potential energy. They parameter- 
ized the mechanical production rate in terms of the friction 
velocity but neglected the viscous rate of dissipation and 
the effect of entrainment shear production. From the 
examination of observations. Turner (1969) deduced that a 
substantial fraction of the part of the work done by the 
wind which goes into the drift current is eventually used 
to deepen the surface layer. 

Miropol ' skiy (1970) and Denman (1973) assumed that 
dissipation is a fixed fraction of mechanical production 
and the remaining turbulent kinetic energy goes to downward 
buoyancy flux. Pollard, Rhines and Thompson (1973) apply 
the slab approach to the oceanic mixed layer but they com- 
plete the entrainment problem with a different mechanical 
energy requirement. In this model, the energy for the 
entrainment is derived directly from the mean flow and the 
intensity of the turbulence is not considered to have an 
active role in the entrainment. Niiler (197^) in his three 
layer model which is a combination of Kraus and Turner, and 
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Pollard et al considered that both turbulent kinetic 
energy and the mean kinetic energy are important in the 
entrainment mechanism. Although he added the entrainment 
production term, the surface production minus dissipation 
was still parameterized in terms of u* . Therefore, 
dissipation was not permitted to adjust to other varying 
atmospheric conditions or to entrainment shear production. 
A more realistic parameterization for dissipation is not 
possible without actually computing the mixed layer turbu- 
lent kinetic energy intensity <£> . It cannot be done with 
only the surface scale u^j.. 



D. REVIEW OF THE RECENT LITERATURE 



1 . The Garwood Model 

Garwood (1976, 1977) developed an ocean mixed layer 
model using the Navier-Stokes equation of motion with the 
geostrophic component eliminated, the continuity equation 
in incompressible water, the heat equation from the first 
law of thermodynamics, the conservation of salt equation, 
and a linearized equation of state. A summary of the 
modeled equations follows; 

entrainment buoyancy flux. 



'laW't-h) 



< E> 
h 



( 1 - 1 ) 
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budget for horizontal components of turbulent 
kinetic energy, 



^ -(<£>- 3 <^» < E>" 

- f C <!>''■ 4 

budget for vertical component of turbulent kinetic 
energy, 

4 4 (<s>- 3 <^»<E>“ 

2 dt ^ 2. ^ 

- iC<E>"+^h)<S> <^‘3 

conservation of mean buoyancy and mean momentum 



hi^r-t/W'(o) + b'vNy'(-h) 4 adz 

d t y-h 

i-(h<c>)= - c'w' ( 0 ) -i-f<2>h 

dt 

"jump conditions" at bottom of mixed layer 



( 1 - 4 ) 

( 1 - 5 ) 



,-Fw'C-h)= AC gf A 



(1-6) 

( 1 - 7 ) 



with the following notation: 



a- 3 Ui 4 U;' 

C = U 41V 
H a U'* 4 V'’4- 



Lija I C'W'Co^l 



I ^ (^-C) 

X 






U< 5 |* 



( ) = 



1 | 7 W 

l.-»a 9 



^ Va r V2 
Xi/2 j-Va 



( 






' I : ah/at X o 

■ 0 : ah/a+ <0 



dv dy 
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To derive the equation (1-1) , the equation of the 
turbulent kinetic energy per unit mass was derived from 
the equation of motion using the Boussinesque approximation, 



at 



( I) It, (I'^ o - It, (“Tl ) + i, urjj 



3Xj 






aX; 



. Fauj , ^ a U J. nrp r. 

+P ^ +Urb<i,-: 



(1-8) 



e. ax; ^x;^ ■ sxj/ ' sxj aXj 

Using the boundary layer approximation and the continuity 
equation, and assuming that viscous diffusion and dissi- 
pation of the mean kinetic energy are negligible, the 
equation (1-8) was written: 



It = +v'w'fj+-gw^-e 



where 



(1-9) 



P - -0 — ( U' Siii ) 



With the assumption that the turbulence of the overlying 
mixed layer provides the energy needed to destabilize and 
erode the underlying stable water mass. The local (at z=-h) 
turbulent kinetic energy budget is the basis of the entrain- 
ment hypothesis. Assuming a quasi-steady state of z=-h, 
and that within an active entrainment zone, the most 
significant source of energy for mixing is the convergence 
of flux of turbulent kinetic energy — ~ Fw' ( j|. +^)] , equation 
(1-9) was reduced to: 






( 1 - 10 ) 
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Estimating the time scale, Tk • required to transport some 
of the turbulent kinetic energy ,(E>, to the viscinity of the 
entraining interface, 






- iBi 

-h 



( 1 - 11 ) 



and assuming that the mixed layer depth is proportional to 
distance over which turbulent kinetic energy must be 
transported by the vertical component of turbulent velocity, 
Te was taken to be K/<w^and therefore from (l-lO) 
and (1-11). 

gw'c-h) 

To get the equation (1-2) and (1-3) > the equation (1-9) 
was first divided into the horizontal and vertical components 
by assuming an isotropic dissipative structure and recalling 



pt / •> 

that ^ — =0 IS included in the equation (1-9) » 



1 9(u'^4y**) 

2 

and 



'r- 2 . _ag a-12) 

3Z L 2. J e, 3 



(1-13) 



1^. 9. r ,(«:> + £;)] +£'8W' +-Ev_Xfe 

a -at 52 . L ' z e, 'J ^ ft. 3 

the equation (1-12) was integrated over the range from r=-h-f 
to "2 -o (Figure 1) 



1 r 0 9(V*-+V‘*i 


dz r-C® , 


^ 1 




2 9t 




9Z 1 


, V 2 . /J 






(1-14) 
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Applying Leibnitz's rule and using the assumption 



r o a ( 

3t 






(1-15) 



Defining a dissipation time scale as =11^ under 

the assumption of locally isotropic turbulence and con- 
sidering that 4< Z^x in which Xfci is the 

convective time scale of large eddies (proportional toh/<E^ 
and is the time scale of planetary rotation 

(proportion to f** ), therefore c X and 



r° feJz =<£>*'*+ 



(1-16) 



The vertical integration of the pressure redistribution term 
is an important source or sink term for the individual 
turbulent kinetic energy budget even though ^ ^ =0. 

Following Rotta (1951) i Lumley and Khajeh Nouri (197^) » the 
equivalent bulk formulation was 



r ® Pj ^U; 

Therefore , 



iz :zE>* (‘i> -3<“t’>) 



(1-17) 

(1-18) 



Following Kraus and Turner (1967) » Denman (1973) and Niiler 
(1975) > "the shear production term was parameterized: 

(1-19) 

\ 

wherelJ*^! "excess" surface mean velocity in the 

direction of the wind stress. In this case, the inhomogeneity 
of the mean velocity field cannot be neglected. 
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Ij,-? [ M 






is proportional tou*i then (1-18) may be combined 
with the parameterized net input from breaking waves, less 
the loss to radiating internal waves; 



-Cf ( 4 TO' 1+ V-V' g 1 ^Z=^U.^ +'4'“^ (1-20) 



From (1-15) » (l-l6) , (1-18) and (l-20) , (l-l4) was written as; 

L ^ *m u» +1^ ^ C< H > -3 < - 1 (< f h )< 1> d - 2 ) 



Similarly integrating (1-13) over the range from z.c4i-r to r=o. 






+ w'n e<iz 

1-h-iT ^ 

and assuming w'=o at the surface and ^^-^1 negligible 

6# > -h-0 

the vertical component of turbulent kinetic energy budget 
was parameterized by using (l-l6) and (1-1?); 



i («)+ yw'(-h)]4<£>"«l)-3<5P>) (1-3) 

“3(<e/*+ fh)<i> 

From the heat equation, conservation of salt equation 
and equation of state, the buoyancy equation was formed. 
Averaging this buoyancy equation and using the boundary 
layer approximation, the mean buoyancy equation was obtained: 

+ <±1 0 

^ X ft (' 1 OT ^ 
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To get the equation (1-4), (1-21) was first integrated 
over the range from to 2=o 



Using the Leibnitz's rule and taking the limit -^-yo 
the term on the left hand side can be written; 

= |£(h<t»-'&(-h-n II 



( 1 - 22 ) 



(1-23) 



Assuming b^w'(-h-^)=0 using /lb- <b> -b(-h-^) , (l-2l) was 
written; 






(1-24) 



To get the jump condition (1-6) , the equation (1-21) was 
integrated over the range from 2=-h to zaO and from z to 
2«0 « and these two integrated equations were solved 

. Then the result was 



simultaneously to eliminate < 4> 

(z) = pw'(o) 0+ -(Fw'(-h) 






(1-25) 



Similarly the integration of the equation (1-21 ) over the 
range from rs-h-i* to z=-h with dropping the negligible solar 

radiation term and assuming Vw' (-h-5* ) = 0 and approaches 

h 

zero led to 



9idz = -b'w'(-h) 



(1-26) 
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By applying the Leibnitz's rule, 



CM = rt - C ^ 

and by assuming l‘^<lcf<b» — n and Ir*"^ 4^ = ^ 

S i Jz= [b(-b^-lC-h-S)] ^ (1-28) 



Therefore from ( 1 - 26 ) and (1-28), the jump condition (1-6) 
was obtained 




and using this relation the equation (1-24) can be written 
in the form of the equation (1-4) . 

Using the boundary layer approximation and dropping 
the negligible viscous terms, the equations of mean 
motion were combined into a complex form by using 2:UflV 



ac 

at 







(1-29) 



To get the equation (1-5). (1-29) was integrated vertically 
over the range from zs-h-f to z-o with the assumption of 

(-h-^) =0 ? 

|^(h<c>)z-iih<c> Co) ( 1 - 5 ) 

The jump condition (1-7) was obtained in a similar way as 

( 1 - 6 ). 
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The most important aspect of this model was the 
ability to predict the year-around evolution of the mixed 
layer depth together with bulk properties. This model 
provides not only structural prediction but also predicts 
a changing layer depth because the entrainment model pro- 
vides boundary conditions at the moving density interface. 

In addition to the non-linear dependency upon stability, 
he found a Rossby number - - dependence for the 

entrainment rate in his nondimens ional solution. This 
makes possible a cyclical steady state for the boundary 
layer without requiring unrealistic values of upwelling or 
lateral advection in a long term integration. 

On short time scales of the order of a few days, 
the upper ocean, even at higher latitudes, exhibits 
significant barocUnc activity. The importance of shorter- 
period fluctuations in the surface buoyancy flux in 
modulating the long-term response has demonstrated the 
need to know the typical daily heating/cooling cycle, solar 
radiation, evaporation, conduction and back radiation as a 
function of season and geographical coordinates. Some of 
the most fruitful applications of this model lie in the 
link-up with other models, for example, joining with an 
ocean circulation model. 

In summary Garwood suggested in his model that first, 
the planetary rotation is assumed to influence the dissi- 
pation for deeper mixed layers and enables a cyclical 
steady state on an annual basis. Secondly, the rate of 
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entrainment for the stable regime is not accurately 
reflected by a linear extrapolation of the unstable 
situations. This is particularly important in modeling 
the ocean boundary layer. Unlike the atmospheric 
boundary case, most of the solar radiation does not 
penetrate the layer. Therefore, downward turbulent heat 
flux in the oceanic boundary layer is as important as the 
upward flux during the course of both diurnal and annual 
cycles. The non-linearity of interface entrainment 
tendency parameter, which is greatest for stable surface 
boundary conditions, results in a modulation of the long- 
term trend of mixed layer depth by the diurnal component 
of surface heat flux. Thirdly, in this model buoyant 
production is somewhat more efficient than shear production 
as a source of energy for vertical mixing because of its 
unique effect on the vertical component of the turbulent 
velocity. 
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2. The de Szoeke and Rhines Model 



The de Szoeke and Rhines (19?6) study of initial 
mixed layer deepening is reviewed here because the results 
will be critically analyzed in a comparison with the 
Garwood model. The solution to impulsive initial deepening 
is most relevant to this research on high frequency response 
of the upper ocean. Following Niiler (1975) » de Szoeke and 
Rhines (1976) investigated the turbulent energy balance of 
a wind-mixed layer (in the absence of solar heating) by 
adding a new term which represents the energy to spin-up 
the level of turbulence in the increment of the mixed 
layer dh in time dt. The equation is 

IA> CB) (i>> 

where term (A) represents the entrainment damping, term (B) 
represents the shear production at the base of mixed layer, 
term (G) represents the energy needed to spin-up the level of 
turbulence and term (D) is the turbulence production minus 
dissipation, in which it is assumed that the dissipation is 
a fixed fraction of turbulent production. 

The results (using constants and parameter Tn*=Co=l > 
U*s| , -fslcr^sec"* , ps s |o'^ ) were 

i) For slightly later than the initial instant, since 
CfisH ^ 1-^* , the shear production, term (B) , can be 

reduced to u^^^t^ and is negligible compared with other terms. 



28 



I 



I 

I 

I 



Assuming a negligible entrainment damping at this early 
stage, the balance (C)'^(D) holds for small t. Initially 
it is a linear deepening, i.e. , ^ oc u*. and this lasts 
for a fraction of a Brunt-Vaisala period (-^100 sec). It 
attains a depth of 2'^?m while the level of turbulence in 
the incipient layer spins up. The term (C) here represents 
a significant demand on the available turbulent energy only 
during the first Brunt-Vaisala period after initiating the 
mixing. After that, it may be safely neglected. 

ii) For large t (about 100 sec to 1 hr) , a balance holds 
between (A) and (D) . After it starts deepening by turbulent 
energy, the entrainment damping term grows. A t'^’ deepening 
emerges as the turbulent work done by the wind at the 
surface erodes heavy fluid from below and mixes it 
uniformly through the surface layer, raising the potential 
energy. This one-third power deepening proceeds from the 
same mechanism that was suggested by Kraus and Turner 

(1967) 

iii) For still larger t (about 1 hr to 12 hr) , a balance 
(A)~(B) holds approximately. The mean flow has been 
accelerating until after 20 to 120 mins at a depth of 5' — ^20m, 
the production of turbulent energy by the shear at the base 

of the mixed layer dominate other sources and a t'^* deepening 
takes over. This h«t*^* was proposed by Pollard, Rhines and 
Thompson (1972) . 

iv) For later times, the Coriolis force swings the mean 
motion vector away from the surface stress direction. 
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reducing the entrainment shear production (term B) below 
4 2 

the rate u* t assumed above. As t increases substantially 
beyond t= TX/i (corresponding to the lapse of one-half 
pendulum day) . An approximate balance between (A) and (D) 
is established again, leading to t*^^ development. 
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II. GENERALIZATIONS OF ENTRAINMENT BUOYANCY 
FLUX EQUATION OF THE GARWOOD MODEL 



Garwood (1976, 1977) developed the entrainment equation 
from the total turbulent kinetic energy budget by assuming 
that the turbulence of the overlying mixed layer provides 
the energy needed to destabilize and erode the underlying 
stable water mass. The total turbulent kinetic energy 
equation (1-9) at z=-h is: 



t If 1-K = - Ih 



(2-1) 



The importance of the unsteady term in the atmospheric 
energy budget was demonstrated by Zilitinkevich (1975) in 
the case of weak stratification, and also suggested by 
Garwood (1976) for the oceanic planetary boundary layer. A 
more detailed study of the importance of this term will be 
in the next chapter. For the order of magnitude of the 
shear production term, Garwood (1977) found that shear 
production is a fixed fraction of buoyant damping in the 
entrainment zone. Zeman and Tennekes (1977) also suggested 
in their atmospheric planetary boundary layer model that shear 
production is absent if the wind shear across the inversion 
base is zero. Even if there is substantial shear, the 
effect of this term appears likely to be relatively small 
(Mahrt and Lenschow, 1976) , because mechanical production 
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of turbulence tends to be largely balanced locally by 
viscous dissipation through the energy cascade mechanism 
(Tennekes and Lumley, 1972). Since the position ( z=-h ) 
of the entrainment zone of the oceanic planetary boundary 
layer is defined by the location where the turbulent 
kinetic energy vanishes, local dissipation is also assumed 
to be a fixed fraction of the flux convergence term, should 
it be significant in the energy budget. Therefore, follow- 
ing Garwood (1977) » a balance holds among the flux 
convergence term, buoyancy damping term and unsteady term. 
The equation (2-1) is reduced to: 



-L §1 

2. 




(2-2) 



Since the boundary erodes and retreats depending on the 
level of the intensity of the turbulent kinetic energy, 
it may be assumed that the level of turbulent kinetic 
energy at the moving interface remains unchanged in the first 
approximation. Therefore, 

gjj I 

<it ‘“H”" Bt St Lh 

or (2-3) 

3^1 - - Hi I 

ilEl-h az at l-h 

Assuming that the turbulent kinetic energy decreases 
linearly with depth in the mixed layer and using the finite 
difference approximation, (2-3) can be written: 
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-h 



- - fl. ^ 2^> Sh 



^z. at 



-^-(-h) at 



h at 



A 



(2-4) 



where <E> is the average value through the mixed layer. 
Employing the parameterization (l-ll) , 7!e = ‘*'he 

jump condition (1-6) , 

j, z ~ ^ Fw^ (-h) s 

z n at h A & h 



and this forms 



(-h) z (2-5) 

h +<E>/a1 

This more general equation of entrainment buoyancy flux 
replaces the equation (1-1) of the section I.D.l for the 
final closure. 

Therefore, the entrainment buoyancy flux, 
and the rate of layer deepening from (1-6) are functions 
of the bulk values, or overall layer averages, of both 
the total turbulent kinetic energy and the relative dis- 
tribution of this energy between the horizontal and vertical 
components. The advantages are two-fold: 

i) The rate of deepening is an integral property--not 
just a function of the local properties of the turbulence 
field, which may be vanishingly small at 2s-h inspite of 
a significant rate of deepening. 
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ii) The relative distribution of turbulent energy among 
the three geometric components is taken into account. 

Being able to include the root-mean-square vertical 
velocity, is particularly important in realistically 

treating layer "retreat", ^ <o. 

dt 

In the case where h aL »<E> , equation ( 2-5) 

returns to (1-1) of the original model. On the other hand, 
in the case < E> » which is the case of weak 

stratification comparing to turbulent intensity (2-5) 

sh 

reduces to s . That is, if the stratification 

is close to the neutral lapse rate, the entrainment proceeds 
at a rate proportional to the root mean square vertical 
velocity, < . The comparable situation in the Zeman 

and Tennekes model is that 6'w' tends to zero as the lapse 
rate, 7” , tends to 0 (adiabatic conditions). Thus the 

energy budget in the absence of shear term and dissipation 
terms becomes a balance between the unsteady term and 
flux. .convergence term i.e., — oc CZ where is 

the average standard deviation of vertical velocity within 
the mixed layer and will scale with < 
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III. DISSIPATION AND FREQUENCY RESPONSE 



A. DISSIPATION AND ITS TIME SCALE 

Turbulence is always dissipative and this dissipation 
occurs through viscous shear stress doing deformation work, 
which increases the internal energy of the fluid at the 
expense of kinetic energy of the turbulence. When these 
local losses due to dissipation are balanced by the local 
inputs of energy, then the turbulent structure in a given 
shear flow might be in a state of dynamical equilibrium. 

Since the turbulence consists of fairly large velocity 
fluctuations governed by nonlinear equations, one may expect 
a behavior like that exhibited by simple nonlinear systems 
with limit cycles. In a rapid energy transfer situations in 
which past events do not dominate the dynamics, one may 
expect that this limit cycle type of equilibrium is governed 
by the local length scales and time scales of the mean flow. 

In this turbulent motion the largest eddies may be as big as 
the width of the flow, which may then be an appropriate length 
scale in the analysis of the interaction of the turbulent and 
the mean flows. Large eddies lose a significant fraction of 
their kinetic energy within one "turnover time", but a 
negligible fraction is directly dissipated. This implies that 
the nonlinear mechanism that makes small eddies out of large 
ones is as dissipative as its characteristic time permits. 
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The nonlinear mechanism is ultimately dissipative hecause 
it creates smaller and smaller eddies until the eddy sizes 
become so small that viscous dissipation of their kinetic 
energy is dominant (eddy Reynolds number of order one). 

Since small-scale motions tend to have small time 
scales, one may assume that these motions are statistically 
independent of the relatively slow, large scale turbulence 
for a sufficiently large Reynolds number for the mean flow. 
Therefore, the small-scale motion should depend only on 
the rate at which the large-scale motion supplies energy 
and on the kinematic viscosity. Thus it is assumed that 
this rate of energy supply should be equal to the rate of 
dissipation. If the dissipation rate, € , can be related 
to the length and velocity scales of the large scale tur- 
bulence, one can assume that the rate at which large eddies 
supply energy to small eddies is proportional to the 
reciprocal of the time scale of the large eddies and 
therefore viscous dissipation of energy can be estimated 
from the large-scale dynamics which do not involve viscosity. 

The amount of kinetic energy per unit mass in the large- 

scale turbulence is proportional to and the eddy tim.e 

scale is assumed to be proportional to where 1 

represents the size of the largest eddies or width of the 

flow. From this general idea, Garwood (1976,1977) derived 

, ^ 

a total dissipation time scale, s -V (for 

n 
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derivation see section I.D.l and for more detail see 
Garwood (1977)) ty including the planetary scale of 
dissipation. This is the integral time scale of the 
turbulence. Depending upon the time scale ( Tf ) of the 
surface forcing (variation in atmospheric conditions), an 
interaction between two time scales might be expected. 

B. SUBJECTIVE EXPECTATION AND ANALYTICAL SOLUTION TO 
INTERACTION BETlffiEN Tk AND Zp 

1 . Constant Dissipation Time Scale 

Garwood (1976) developed a solution for forcing 

period much longer than the integral time scale of the 

turbulence. This solution is governed by the dissipation 

time scale ( ) . For these time scales one can assume a 

quasi-steady state so that (1-2) and (1-3) become diagnostic 

equations for < w'*> and < u'* +v'^> . In this way, the 

entrainment buoyancy flux (Fig. 3) can be computed directly 

as a function of the surface fluxes of momentum and buoyancy. 

2 

However, fluctuations in surface boundary conditions u^ 
and b'W'{6) of sufficiently short period require the con- 
sideration of the unsteady terms . 

A very simple model which can be solved analytically 
provides an indication of the importance of this unsteadiness 
and how it might be treated. In the simple model, only two 
large terms, shear production and dissipation, are balanced 
by the time rate change of the turbulent kinetic energy. 
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(3-1) 



+ iMl - yi. ( I + ^ Sir\ tut ) 

St Tt h 



where "a" is the amplitude of fluctuation and U)= 2Ti/rp 
as its frequency. Assuming that the dissipation time 
scale is constant and defining a dimensionless time, r=t/T^ 
a dimensionless turbulent kinetic energy, (|J = 




f 



t 



and frequency, ou = 2iT T^/ tp 



( 3 -I) is written; 



^ + (b = I ‘as\'n(a;»t) 

at ^ 



( 3 - 2 ) 



Eor cyclical steady state, the analytical solution to (3-2) 
is 






( Sih (ju^r — Cos ) 



(3-3) 



For the case a=o , =1 which is the non-cyclical 

steady state solution. For finite amplitude, a, the full 
solution ( 3 - 3 ) is needed. In terms of frequency, for the 
low frequency fluctuations, (i.e., long period 

compared with ) , the approximation of quasi-steady 
state is satisfied. Depending upon the Rossby Number 
(Rq = *l»/(fh)) and the mixed layer stability parameter 

(^* ~ iT^ ^ ' this low frequency fluctuation corresponds 

to a fluctuation period of several hours or longer. From 
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the complete solution (3“3),‘the important aspect is the 
relative magnitude of the response associated with the 
fluctuating component (fig. 4). For very high frequencies, 
» 1. the amplitude of the response is negligible. In 
other words, the high frequencies are filtered out. This 
becomes a possible cause of error in using observed winds 
to drive a model having an integration time step smaller 
than or comparable with • If the quasi-steady state 
assumption is made to facilitate obtaining a solution, but 
without filtering or smoothing the surface boundary 
conditions, an incorrect high-frequency response will not 
only be present but it may bias the mean trend. 

2. Non-constant Dissipation Time Scale 

In this case the same simplified equation (3-1) is 
considered, but the dissipation time scale is no longer 
constrained to be constant. Using the convective dissi- 
pation time scale 7^= ■ of Garwood (19?6)^ (3-1) becomes 



e<E> _ 





( 1+ a SiTJ Wt) 



(3-^) 



To examine the significance of the unsteady term, equation 
(3-^) is first nondimensionalized by defining 

T , 



and 



= 



(juh 

u* 
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then (3"^) is transformed to: 



I +a Sin ( CtWT) 
9T 



(3-5) 



The equation (3~5) is solved numerically for several 
combinations at amplitude (a=0.1, 0.5 and 1), and angular 
frequency of forcing ( - 200n , 20 tl , 2 IT ,0.2TI, and 

0.02TT). Initial adjustment was achieved at about . 

There is also a tendency that the higher the forcing 
frequency , the shorter the initial adjustment period 

that is needed. At the early stage of the wind forcing, the 
change of turbulent kinetic energy is mostly due to the 
wind forcing itself. As the energy supply from the wind 
continues, and after an initial adjustment is achieved, the 
balance among the three terms is changed and becomes a 
function of the amplitude and frequency of the forcing. 

Different responses of the unsteady term 
to different frequencies of forcing are shown in figure 
5 a-b for the case of small (a=0.1) and large (a=l) 
amplitude fluctuations. For both amplitudes, there is a 

tendency for the mean value of the unsteady term, , 

L4* 9t 

to vary with frequency even though it is very small. This 
is seen more clearly in the case of small amplitude fluctu=- 
ation (a=.l) because of the magnified scale for the small 
amplitude case. The mean value is a little above zero for 
high frequency and it approaches zero as the frequency is 
getting lower. This effect is also shown in figure 6 which 
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are plots of the case 1 ); Sincot function of the wind 
forcing versus lUt , case 2) ; fluctuating component of 
dissipation normalized on the amplitude of the wind 
production ( ) versus u»t for a=l , and 

case 3); same as case 2) except a=.l. Comparing the 
fluctuating component of dissipation, case 2) and 3) » with 
the constant wind forcing case l) , for each of three fre- 
quencies, one can see that the higher the frequency, the 
smaller the fraction of wind production that is dissipated 
in the mean. In other words, more energy is stored in the 
mean to deepen the mixed layer. 

Figure 7 shows the peak- to -peak amplitude of the 
unsteady term " 9^ versus frequency for three 

amplitudes of fluctuation a=l, 0 . 5 aJ^d 0.1. For each case, the 
peak-to-peak amplitude is almost a constant close to 2,1 
and .2 respectively over the frequency range - ztt to 

Zoo IT is negligibly small for the frequency 

range below -OZTT . Some interesting transition 

cases occur over the frequency range ^ 02TT "to 2 tt 

in which the adjustment time scale and forcing time scale 
are comparable, and therefore interact. For relatively 
low-frequency surface forcing, ^ .q^TT (nondimensional- 

ized period Tp = ^ ^ 100) , regardless of magnitude of 

oscillation amplitude, the turbulent kinetic energy budget 
is in a quasi-steady state, i.e., dissipation is nearly 
proportional to the wind shear production as indicated by 
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the negligible peak-to-peak amplitude of unsteady term, 

. Therefore, for this case, ( 3 ~^) reduces to 

at 

<i> = u» ( j Sin (3-6) 



For relatively small amplitude oscillations (a<, 0 . 5 )» 
regardless of frequency, it may be reasonable to assume 
a quasi-steady state for many cases in which a high degree 
of accuracy is not required. But for relatively large 
amplitude oscillations, the unsteady term is more clearly 
frequency-dependent. For a < 0.5 and in the transition 
frequency range below c.2iT > it may not still be 

a bad approximation to assume a quasi-steady state and use 
(3-6) in an approximate sense. On the other hand, for 
a i . 5 and frequency = . 2iT , this unsteady term is 

U4 

no longer negligible. Therefore it is strongly suggested 
that the full form of equation (3-^) be used to avoid an 
incorrect high frequency response. 

Figure 8 shows the phase difference between the 
wind forcing and the dissipation as a function of the 
frequency of the forcing. For the frequency range 
there is almost no phase shift. As shown by the negligible 
amount of hysteresis in figure 9a, the dissipation is nearly 
proportional to the wind shear production after the initial 
adjustment period o< <4 .In the other extreme of 
high frequency forcing ( ^ > 2.0 tt ), there is almost a 90° 
phase shift between wind production and dissipation. As 
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shown in figure 9d, the adjustment period is longer 
relative to the forcing time scale, but dissipation 
converges to a constant value that is independent of the 
frequency and magnitude of the forcing oscillations. For 



time scales are comparable, the phase shift increases 



9b and 9c show that since the adjustment time scale and 
the forcing time scale are comparable and interact, this 
causes the ellipsoid to achieve a maximum width in this 
domain, 1< Tp< 10. An important result is revealed by 
figure 9"b which shows that the amplitude of the response 
cycle is nearly proportional to the amplitude of the 
forcing cycle, inspite of the hysteresis. Therefore, if 
phase information is not crucial, it is accurate for this 
case to assume a quasi-steady state. 



the frequency 




, in which the two 



nearly from zero at .02TT to nearly 90° at 20Tt . Figure 
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IV. ASYMPTOTIC PROPERTIES OF THE GARWOOD MODEL 



A. NUMERICAL RESULTS OF GARWOOD MODEL 

If quasi-steady state is not assumed, the full prognostic 
equations may still be solved numerically. The full 
equations of Garwood model, (2-5)» (1-2), (1-3) » (1-^) and (1-5) 
with the jump conditions (1-6) and (1-7) are solved 
numerically for about a one-day period using hypothetical 
initial conditions, h=lm, <e> = • o 3 cmVsec* , < w’*> = . ol , 

Q:o , v~o and three different cases of , 

5xio"^ and )• size of indicates the 

degree of stratification and is related to the lapse rate , T , 
as follows; 

?W'(oJ <Jt (4-1) 

Z n -'o 

where P = ^ C: s *25 ^ for the simple case of 

^21 3Z 

vertically homogeneous salinity. Assuming the surface 
buoyancy flux to be negligible, » 4b = BXioT^ Cm/sec^ 

corresponds to = I *c/;t 5o rn ( i . e . p = \ o'® sec~^ ). 5Xlo^Cw/se<f 
corresponding to (i. e. p- »o“^ sec'^ ). The 

constant employed here is . To simplify the problem, 

the momentum flux u'w'(o) - 1 cmVscc’, v'w'so , the constant 
frictional velocity | cm/ sec ,A=o,17w'=o and -^=10”^ see"* 

are also assumed. Using these conditions, three different 
curves of each parameter are plotted for each lapse rate. 
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Figure 10 shows that the stronger the stratification, the 
larger the amount of turbulent kinetic energy in the mixed 
layer. It starts from ^ .03cmV^*and very rapidly 

increases during a short period immediately following an 
impulsively applied surface forces. 4E > shows a maximum 
at t=3~5 hrs. , a minimum at t=17— ' 18 hrs. and a second 
maximum at t=24 — ^25 hrs. The figure also shows that the 
weaker the stratification, the earlier the maximum is 
reached. The vertical component of turbulent kinetic 
energy shows a similar pattern except for the case of 
weak stratification ( ^ = l *c/aSo wj ) in which it decreases 
monotonically (fig. 10b). Figure 11 shows versus 

time for and • Starting with 

0.33f ii reaches 0.15» 0.14 and 0.13 respectively 
after a one-day period. These sudden reductions in <w' ^>/£> 
occur at a very early stage and then stay nearly constant. 
In other words, vertical components of turbulent kinetic 
energy relative to total kinetic energy adjust rapidly in 
the early stages, and, the stronger the stratification, the 
larger will be the vertical component relative to total 
kinetic energy. 

Figure 12 shows changes in At with time. Each 
increases monotonically, although more rapidly in the case 
of stronger lapse rate. Larger values of Ab contribute to 
more turbulent damping due to mixing at the interface. 
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Figure 13 shows the change in current with time, and 
indicates that the stronger the stability, the larger will 
be its magnitude. Vlith the imposed east wind, the currents 
start flowing westward at the initial instant. The west 
component reaches a maximum value of 17 cm/sec, with the 
northward component current still growing. As the west 
component decreases to zero, the north component current 
reaches a maximum of 25 cm/sec. This slowly turning current 
direction is caused by Coriolis force which deflects the 
current direction to the right in the Northern Hemisphere. 
After one inertal period, both u and v decrease to zero. 

The m.aximum of + occurs at about t=4hrs and the 
minimum is at about t=17hrs . Also the figure I 3 shows 
that the magnitude of each successive maximum value decreases 
with time. 

Figure 14 shows the monotonically increasing mixed layer 
depth with time. All cases start from Im and after lO^sec 
the mixed layer depth is 45m for ^ = 1° /250m, 23 m for 

= l°^/25m and 11m for = l°*^/2.5ni. These indi- 

cate that the larger the stability, the lesser the deepening 
it is. 

Figure 15 a-c show (or 

because UnS. Icny^Bcis assumed) and for each lapse 



rate. The middle curve of each figure indicates that 



dh 

^3t 



decreases more rapidly in large stability than in weak 
stability. For each case, and ^ 
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are not parallel prior to about one inertial period, i.e., 

^ oc and ^ oc would not produce the same deepening 

even though |yw'(®) = o is assumed. The figures imply that in 

general the deepening is not proportional to 

oru'^, but rather depends upon the full form of entrainment 

equation 4^ ~ . An exception occurs if h^toro 

cH: + 

so that during the early stages (Figure 15d) . 

Figure l6 shows the inverse of the bulk Richardson 
number changing with time. For t < 20 mins, RT* is large 
in the case of weak initial stability because of the small 
Ab . After 20 mins, the larger values of occurs for 
the cases with larger initial stability. This is due to the 
strong mean shear developed at this stage. At one inertial 
period, R;”* goes to zero because of the stagnation in 
current at this moment. Halpem (197^) suggested that if 
local Richardson number ( R?i = ^ / ) is less than 

one-fourth, instability may occur, and if R;i =1, marginal 
stability exists. He also suggested in his experimental 
study of Northern Pacific Ocean that is a critical 

condition for the occurence of shear instability. Therefore, 
for the case of ^ r i**c/a 5 o m , there is an instability at 
about t=100sec. If Ru is considered to be an average of 
over the mixed layer depth and therefore , 

there may also be an instability at t=4hrs for the case of 

|«»c/2.5m • This instability leads to a maximum turbulent 

kinetic energy at this time interval. 
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Therefore, how the turbulence develops with time, as 
indicated by equation (1-2) and (1-3) > can be summarized. 
The total turbulent kinetic energy (fig. 10) <£> jumps 
due to impulsively applying surface stress to an initially 
resting fluid, whereas the energy produced by any mechanism 

w 'S/X 

is dissipated at a rate proportional to . , Then, the 

shear production is increased to a maximum value even 
though the damping also increases. The maxim.um turbulent 
kineiic energy exists at about t= 4 hrs, and then decreases 
toward a minimum mainly due to decreasing shear production 
at about one inertial period. 

The more detailed, term-by-term study is left for the 
next section. 

•B. ASYMPTOTIC REGIMES IN MIXED LAYER DEEPENING 

Combining equations (1-2) and (1-3) of the Garwood 
model, the total turbulent kinetic energy budget equation 
is obtained. 

2Yn ui + + HFw'(o) - 2 fh<B» 

Using ' one can write the above equation in the 

following order; 

- 2fn uj -2 4 -?ke» ( 4 - 2 ) 

(A) C6) C.C) CD) 

The corresponding equation from de Szoeke and Rhines (19?6) 

p 

may be obtained by dividing equation (I-30) by h . 
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(4-3) 



^ ^ ^ ^ ^ ■+ Co Ui# ” 2 jn, Ujf 

CA) (.8; CC ) 0> 

Equation (4-3) is originally the Niiler model except for 
the term (C). The term (A) in both (4-2) and (4-3) is an 
entra.inment damping term which is a upward buoyancy flux 
at the base of this mixed layer. The term (B) in both 
equations is the shear production term. The term (C) 
in (4-3) is the energy needed to spin up the level of 
turbulence, , in the increment of mixed layer dh in 

time dt, while the term (C) in (4-2) is the nonstationarity 
term which shows the incremental change of the turbulent 
energy budget d(h<E> ) in time dt. The term (D) in both 
equations is the wind shear production minus dissipation. 
Before analyzing the numerical solutions it might be useful 
to examine the characteristics of both equations and their 
differences . 

i) Both models contain nonlinear equations, i.e., the 
Garwood model is a set of six nonlinear differential 
equations which can be solved implicitly for each parameter, 
while the de Szoeke and Rhines model has one nonlinear 
differential equation. 

ii) Both models have numerical problems when the zero- 
initial condition (h(t=0)=0) are used, although such 
conditions are not likely to occur in nature. 
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iii) The nonstationarity term in Garwood, 

h is not exactly comparable to Co 

of de Szoeke and Rhines. Even with the assumption of no 

surface buoyancy flux,<E/u^ does not approach a constant 

. ^=rc/i5m and g =. | up to 
one inertial period. After this period, <E/u^j. = constant 
holds approximately, although it is still not exact except 
for the case of weak stratification (figure 10a). Even 
when the condition<E/u* = constant holds, there is a 

I — 

difference of between the two models. 

iv) Term (D) in the Niiler (1975) model is 
parameterized in the manner; 




where " <1 " is the depth of the perturbation energy prod- 
uction zone near the surface. This implies that only the 
near surface production is dissipated. In other words, 
dissipation is a fraction of production term near the 
surface. None of the shear production energy at the base 
of mixed layer and the surface buoyancy production energy 
(this term is not considered for this study) goes to 
dissipation. Because of this parameterization, the shear 
production at the base of mixed layer has only an interaction 
with the entrainment damping term. On the other hand, in 
the Garwood model, dissipation is parameterized in a way 
that all the production including the entrainment shear 
production and the surface buoyancy flux is dissipated. 
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Term (D) in the Garwood model is the turbulent production 
over the whole depth minus the total dissipation, i.e., 
the viscous dissipation as well as the planetary scale 
dissipation. This second dissipation term is not included 
in the Niiler model. If the planetary scale dissipation 

is neglected and —2 = a constant is considered, then 

•-V 

term (D) in the Garwood model approaches term (D) in the 
de Szoeke and Rhines equation of the Niiler model within 
a constant factor. 

With these characteristics and differences in mind, two 
cases of the Garwood model and one case of the de Szoeke and 
Rhines model are solved numerically using the same conditions. 

Case 1) The equations (2-5), (1-2), (1-3), (l-4)and (1-5) 
with the jump conditions are solved numerically using the 
same constants and initial conditions as in the case of 

) in Section IV. 1 (i.e., 

A Z 

Parameters and boundary conditions : m = g , U^s.\cm/^ec 

— 1 ^ 0 ini'tis.l cond.i’tions 5 

h - I 00 cw , <E> s ♦ o3 = . 0 1 cmVsec* » ~ 5X cm/sec') 

Case 2) is similar to case 1) except that there is no 
planetary scale dissipation term in (1-2) or (1-3). 

Case 3) Equation (4-3) of de Szoeke and Rhines is solved 
numerically using C,=m.=l , = « cw/sec , sec*< 

and p 2 a- 1 0-4 sec . 

Each solution is displayed from t=50sec because the 
initial values of terms (C) and (D) are large. The three 
cases of term (A) are plotted in Figure 17* Maximum 
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entrainment damping occurs at t=4^5 and a second 

maximum at t=25 hrs for each case. The entrainment damping 
of the de Szoeke and Rhines model is twice as large as that 
of the Garwood model at the time of the first maximum and 
remains larger than that of the Garwood model throughout. 
Case 2) of the Garwood model and Case 3) of de Szoeke and 
Rhines are very close at time equal to one inertial period. 
The entrainment damping of Case 2) (i.e., without planet- 
ary scale dissipation) is larger than that of Case l) of 
Garwood model after t=15 mins and the difference increases 
with time. For the time interval less than 15 mins, the 
planetary scale dissipation term does not make any con- 
tribution to the entrainment damping. The large damping of 
term (A) in the de Szoeke and Rhines equation compared to 
the Garwood model occurs inspite of the large shear prod- 
uction in the Garwood Model (Figure 18). This may he due 
to the fact that none of the shear production of the 
de Szoeke and Rhines equation dissipates. Instead, all 
of this energy goes to deepen the mixed layer, and 
therefore the entrainment damping grows larger to offset 
the larger entraining energy. The reason for larger 
damping in case 2) than in case 1) is because there is less 
dissipation in case 2) (without planetary scale dissipation) 
and thus there is more entraining energy available, against 
which entrainment damping grows. Since this planetary scale 
dissipation increases with time, the difference in the 
entrainment damping between the two cases also increases 
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with time. Since the entraining due to entrainment shear 
production decreases to zero (Figure 18) at one inertial 
period in all cases, the net entraining energy in case 2) 
and case 3) is "the same at this time stage and therefore 
the entrainment damping in case 2) and case 3) is almost 
the same. Whereas, case l) has less entraining energy 
because of planetary scale dissipation. This leads to less 
entrainment damping at this same time stage. 

Figure 18 shows term (B) for the three cases. The 
figure shows that the entrainment shear production of the 
Garwood model is not affected by the planetary scale 
dissipation term and the term (B) of the Garwood model 
is larger than that of the de Szoeke and Rhines model at 
least for one-day period. It also shows that the time of 
maximum shear production agrees with the time of maximum 
entrainment damping as discussed in term (A) . The reason 
for the smaller entrainment shear production of the de Szoeke 
and Rhines model is because there is no dissipation 
attributable to the entrainment shear production, and this 
leads to larger deepening in the de Szoeke and Rhines model 
than in the Garwood model (Figure 19)- Therefore, the term 
(B) of the de Szoeke and Rhines model, which is divided by h , 
shows a smaller value of extrainment shear production. 

Figure 20 shows the value of the (C) terms. Again, the 
planetary scale dissipation term of the Garwood model does 
not affect this term. This term in the Garwood model is 
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larger than that of the de Szoeke and Rhines, although 
the difference between the two models diminishes in time. 

The reason for term (C) (Figure 20) being larger in the case 
of the Garwood model , and possible for the decrease in the 
difference with time is because of the omission of the term, 



d<i> 



, in the de Szoeke and Rhines model. . In addition 



this is also because <E>/U# tends to be large at first and 
only becomes nearly constant after about an inertial period. 

The last and important (D) terms are plotted in figure 
21. The term (D) of de Szoeke and Rhines is equal to 2.0 
and thus appears as a straight line. On the other hand, 
the term (D) for case 1) and case 2) of the Garwood m.odel 
varies with time. After several minutes, the planetary 
scale dissipation term starts to act and because of this 
action, the term (D) of case 1) is always less than that of 
case 2) with larger differences at the maximum and minimum 
values. After one inertial period, case 2) and case 3) are 
close together, because there is no entrainment shear, 
therefore no dissipation of entrainment shear production 
(this is always the case for de Szoeke and Rhines) and no 
dissipation from planetary rotation in both cases. The 
larger differences of the term (D) between case 1) and case 
2) at the time of minimum is because of the larger value 
of due to large entrainment shear production which 

contributes to larger planetary scale dissipation, 

The term (D) of the Garwood model shows a negative value 
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between t=1.3 hrs and t=10.1 hrs. The reason for this 
negative value during this time interval is because the 
term (D) is not a net production term, i.e., this 
term does not include the shear production at the base of 
the mixed layer but does include the dissipation from that 
entrainment shear production. Therefore, this term shows 
a minimum during this time interval whereas term (B) 
(Figure 18) shows its maximum. 

To examine the asymptotic regimes, the solutions for 
case l) , case 2) and case 3) are plotted for separate 
(A), (B) , (C), and (D) terms and also for B/A, C/A and 
d/a following de Szoeke and Rhines (see also Figure 2). 
However, since the results for case l) and case 2) are al- 
most the same, only case l) will be considered. From the 
initial instant to about t=5 mins a balance (C)^(D) holds 
for the Garwood model (Figure 22). 



This is similar to the balance in the de Szoeke and Rhines 
model except that the duration for their model is less than 
100 sec. In the de Szoeke and Rhines model 



a(h<£» 

dt 






(4-4) 




or 




( 4 - 5 ) 
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This (C)^(D) balance for both models implies that the 
initial deepening proceeds due to turbulent motion in the 
mixed layer. 

After this period, approximately from t=5 mins to 
t=l8 mins, there is no simple balance in the Garwood 
solution. Therefore, the full form of the turbulent 
kinetic energy equation is required during this period. 
However, in the de Szoeke and Rhines model, a balance 
between (A) and (D) holds rather clearly from the early 
stage to about 1.2 hrs. For the Garwood model, this is 
the period that the value of nonstationarity is decreasing 
and entrainment shear production is increasing, although 
both the entrainment damping and the production minus 
dissipation are increasing to offset the value of non- 
stationarity and shear production. From t=l8 mins to 
about t=12hrs, a balance (A)/^(B) holds for the Gairwood 
model , 

or 

O . L£5L*~| ('*-6) 

hat 

This is the regime in which the bulk gradient Richardson 
numbers is nearly constant (RiCrl ). For this time domain, 
the value of nonstationarity is decreasing and the strong 
shear production is balanced by the entrainment damping. 

In the de Szoeke and Rhines model, this balance (A).~'(B) 
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holds from t=1.2 hrs to 8.8 hrs , which is somewhat shorter 
period compared to the Garwood model. Finally following 
t=12 hrs, a balance (A)'^(D) holds for the Garwood model 

~ 2m 4- 4h<¥>) (4-7) 

Likewise in the de Szoeke and Rhines model, the same 
balance (A)->-(D) holds although this begins earlier at 
t=8.8 hrs. 
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V. SUMMARY 



The one-dimensional model of the ocean mixed layer of 
Garwood (1976, 1977) has "been further generalized. The 
entrainment buoyancy flux equation in this model is derived 
from the turbulent kinetic energy budget at the base of the 
mixed layer. In earlier models, the unsteady term 
has been considered to be negligible in the derivation of 
the entrainment buoyancy flux equation. From the reduced 
form of the turbulent kinetic energy equation in which a 
balance holds among the shear production, the dissipation 
and the unsteady terms, Garwood (1976) nondimensionalized 
and solved this equation linearly using a constant dissipa- 
tion time scale. He suggested that if the quasi-steady 
state assumption is made to facilitate solution, without 
filtering the surface boundary conditions, an incorrect 
high frequency response will not only be present but will 
bias the mean trend. Following these suggestions, the 
unsteady term has been included in the parameterization of 
entrainment buoyancy flux equation. Even though the entrain- 
ment buoyancy flux equation of the Gairwood model is sufficient 
in most cases, this new generalized equation is necessary 
to study the asymptotic regimes. 

To examine the importance of the unsteady term in the 
kinetic energy equation in more detail, the same reduced form 
as the Garwood (1976) equation is nondimensionalized, and 
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solved numerically using a nonconstant dissipation time 
scale , Solutions for variable amplitudes of 

fluctuation and forcing frequencies were obtained and 
their results are as follows. After an initial adjustment 
period o < 4 • 

i) For any amplitude of fluctuations, the quasi-steady 
state assumption may be good, if Wh/U^^ ^0.02TI. 

ii) In the case where the amplitude of the fluctuation 

is smaller than a half of the mean amplitude or UJh/L4<»2TT 
(i.e.,7p> 10) the assumption of quasi-steady state may be 

possible . 



iii) For the frequency tuK/Ut>.2TT , or the amplitude 
being larger than one-half of the mean amplitude, the 
unsteady term should not be neglected. 

Following the lead of de Szoeke and Rhines (1976) and 
Niiler (1975) » "the relative contributions of entrainment 
damping, entrainment shear production, nonstationarity , 
wind shear production and dissipation were examined. The 
specific differences between two models are; 

i) The Garwood model is computed implicitly for each 
term whereas the de Szoeke and Rhines' model is solved 
explicitly. 

ii) Even with assumption..of = a constant, an 

additional term is in the Garwood model for the non- 

stationarity term (C). 
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iii) Dissipation is parameterized in a different 
manner, that is, Niiler (1975) parameterized the 
dissipation term as a fraction of near-surface production 
term, whereas Garwood (1976, 1977) parameterized the 
dissipation term in relating to net production. The 
numerical solutions of both models indicate that larger 
entrainment damping occurs in the de Szoeke and Rhines 
case inspite of smaller shear production compared to Garwood. 
Also, the solutions indicate that a larger nonstationarity 
term occurs in the Garwood model, and that a constant 
production minus dissipation term occurs in the de Szoeke 
and Rhines model compared to large variation in the 
Garwood model including the negative value of this term. 

Following the method of de Szoeke and Rhines, each 
asymptotic regime is classified for both of the models. 

i) From the initial instant to about 5 mins, a balance 
between the nonstationarity and the production minus 
dissipation holds in the Garwood model, i.e.. 



This same balance holds in de Szoeke and Rhines with 
duration from the initial instant to about less than 100 sec. 

ii) From t=5 mins to t=l8 mins, a balance holds for the 
full form of the turbulent kinetic energy equation, whereas 
in the de Szoeke and Rhines model, a balance between the 



d (h^s>) 
<jt 

(C) 



-w zm u/ -2 E>) 

(D) 
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entrainment damping and the entrainment shear production 
holds; 



(A) (B) 



or 



Ri 



i^cl* 

h 



I 



This balance (A)^(B) holds also in the de Szoeke and Rhines 
model from t=1.2 hrs to t=8.8 hrs. 

iv) Following t=12 hrs, a balance between the entrainment 
damping and the production minus dissipation holds; 



i 2 m u/ - 2 ^ li^>) 

(A) (B) 



whereas in the de Szoeke and Rhines model, the same balance 
(A)/^(D) holds, although this balance occurs a little 
earlier, from t=8.8 hrs. 
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Figure 1. Idealized model for ocean mixed layer ( ) . Mixed layer depth 

is (h);(S) is the thickness of the interface or entrainment zone. 
(Garwood ( ] 976) ) 
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Fignre 2. Numerical solution of the complete mixed 
layer model equation of de Szoeke and 
Rhines shov^ing various asymptotic r^fgimes . 
Solid curve: the ratio D/A of terms in the 
equation; dashed curve: B/A; dotted curve; 
C/A. Parameters U;^=l,m» = l, c^ = l (de Szoeke 
and Rhines (] 976 ) ) 
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'igure 3. General solution to entrainment and turbulent kinetic energy 
equations. (Garwood (] 976 ) ) 



64 









I 




J I 1 L__ 

1 rnorilh I day I how 1 mm 

Wtfrt©’' deep, alow-respondirnj mixed layef 

1 1 ! J 

i month 1 (Joy 1 h<xjf I min 

ig^re 4. Response of the mixed layer to 
fluctuations in the surface 
boundary conditions; filtering 
effect. (Garwood (] 976) ) 
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uh/u* = 200Tr , 2it and .02ir, 



66 




Figure 6. Comparison of following fcliree cases for each non— 
dimensionalired frequency (,a ) uh/u j^--2 00 .t , Cb ). -> n and 
(.C ). . 0 2 TT • 

Case ]) sinojt function of the wind forcing vs. cat. 
Case 2) Fluctuating component of dissipation 
normalized in amplitude of wind production vs. wt for 
for a=) . 

Case 3) Same as case 2.) except for a=.] . 
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Figure 7. Peak to peak amplitude of unsteady term 
vs. nondimensionali zed frequency. 




and dissipation vs. nondimensionalized 



frequency 
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ro 9. Nondinensional dissipation versus prescribed wind shear production 
for four periods of osciallation . 
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Figure 10. Cal Tot£il turbulent kinetic energy and 

Cb) Vertical component of turbulent kinetic 
energy versus time for lapse rate 
1 °C/2 . 5ra, ] °C/25m and ]'>C/250m. 
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Figure 11. <w ' /<E> versus time 
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Figure ]2. Ab vs. time for lapse rate 

] ®C/2 . 5m , ] <»C/2 5ni and ]<»C/25Cim- 
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Figure 13. u and v component of mean velocity 
for lapse rate l'’C/2.5m, 
)*C/25m aiid )‘'C/250ra. 
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Figure ]4. Mixed layer depth vs. time for lapse rate 
] °C/2 . 5m , ] ®C/25m and )^C/250rQ. 
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lapse rate (a) ] "0/2 . lira, (b) ) ‘*C/25zi, (c) ] ‘*C/250ia 
and (d) ] ‘*C/250000m. 
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rate ] “C/2 . 5ni , ] '“C/25m and 
] “C/250ra. 
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term in three cases 
Case 1) Garwood — full equations. 

Case 2) Garwood — without planetary scale dissipation. 
Case 3) de Szoeke and Rhines. 




Figure ]8. Comparison of the entrainment shear 
production term in the same three 
cases as in figure ]7. 
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Figure ]9. Mixed layer depth versus time. 
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Figure 20. Comparison of the nons tationarity 
term in the same three cases as in 

figure J7. 




Figure 2]. Comparison of the production minus 
dissipation term in the same three 

cases as in figure 37 
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Figure 22. Numerical solution for the Garwood model 

( ^) in terms of B/A,C/A and D/A. Parameters 
and boundary conditions: m=6 , f = ] 0“ ** sec” ^ , T=N^ 
= ]0“^sec”^ ,u'v' = ]cm^/sec^ ,v''w'-0,b'w" = 0,u^ = ] 
.cm/sec , Ini tial conditions: h=) 00cm , <E>= . 03 
cm^/sec^ , <w ' ^ > = ,0Jcm^/sec^,Ab=5x']0” ^cm/sec^ . 
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Figure 23. Numerical solution for the Garwood model 
without the planetary scale dissipation 
term (a) in terms of A,B,C and D separately, 
and (b) in terms of B/A,C/A and D/A. 
Parameters and conditions are the same as 
in figure 22 . 
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Figure 24 



t (sec) 

Numerical solution for the de Szoeke and 
Rhines' model (a) in terms of A,B,C and D 
separately, and (b) in terms of B/A, C/A and 
D/a. Parameters: m, =c „ = ) , u* = ] cra/sec , f = ] O"** 

Initial condition: 



sec ' , r = N = 1 0" ^ sec , 



. - i r = M 2 = 

h=l 00cm. 
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